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Abstract 

We  observe  that,  under  very  mild  conditions,  an  n-dimensional  space  of  functions 
(with  a  finite  n)  admits  numerically  stable  n-point  interpolation  and  integration  for¬ 
mulae.  The  proof  relies  entirely  on  linear  algebra,  and  is  virtually  independent  of  the 
domain  and  of  the  functions  to  be  interpolated. 

Approximation  of  functions  and  construction  of  quadrature  formulae  constitute  an  ex¬ 
tremely  well-developed  area  of  numerical  analysis;  in  most  situations  one  is  likely  to  encounter 
in  practice,  existing  tools  are  satisfactory.  Much  of  the  research  concentrates  on  obtaining 
powerful  results  under  strong  assumptions  —  designing  interpolation  and  quadrature  for¬ 
mulae  for  smooth  functions  on  subspaces  of  manifolds,  etc.  In  this  note,  we  make  a 
very  general  observation  that,  given  a  hnite  set  of  bounded  functions  /i,  A,  . . . ,  fn-i,  fn 
(either  real-  or  complex-valued)  dehned  on  a  set  S,  there  exists  an  interpolation  formula 
that  is  exact  on  all  linear  combinations  of  /i,  A,  . . . ,  fn-i,  fn,  is  numerically  stable,  and  is 
based  on  n  nodes  in  S  (to  be  denoted  xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn)-  If,  in  addition.  S'  is  a  measure 
space,  and  the  functions  /i,  A,  . . . ,  fn-i,  fn  am  integrable,  then  there  exists  a  quadrature 
formula  based  on  the  nodes  Xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn  that  is  exact  on  all  functions  /i,  A,  . . . , 
fn-i,  fn,  and  is  also  numerically  stable.  Both  of  these  statements  are  purely  linear-algebraic 
in  nature,  and  do  not  depend  on  the  detailed  properties  of  S,  or  of  the  functions  /i,  A,  . . . , 

fn—l,  fn- 

It  should  be  pointed  out  that  all  of  the  statements  in  this  note  follow  easily  from  the 
analysis  found  in  [3] ;  moreover.  Theorem  2  can  be  found  (in  a  slightly  different  form)  in  [4] . 
While  we  cannot  cite  any  earlier  works  where  these  observations  are  published,  it  seems 
unlikely  that  they  had  not  been  made  a  long  time  ago  (perhaps  in  contexts  other  than 
numerical  analysis). 

Throughout  this  note,  S  denotes  an  arbitrary  set,  n  denotes  a  positive  integer,  and  /i,  A, 
. . . ,  fn-i,  fn  denote  bounded  complex- valued  functions  on  S  (all  results  of  this  note  also 
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apply  in  the  real-valued  case,  provided  that  the  word  “complex”  is  replaced  with  “real”  ev¬ 
erywhere).  For  any  n  points  Xi,  X2,  •  •  • ,  Xn-i,  Xn  in  S',  we  dehne  A  =  A{xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn) 
to  be  the  n  x  n  matrix  dehned  via  the  formula 


^j,k  fj^Ak) 


(1) 


with  j,k  =  1,  2,  . . . ,  n  —  1,  n;  we  dehne  the  function  gk  on  S  to  be  the  ratio  of  the 
determinant  of  A[x\,  X2,  ■  ■  ■ ,  Xk-2,  Xk-i,  x,  Xk+i,  Xk+2,  •  •  • ,  Xn-i,  Xn)  to  the  determinant 
of  A[xi,  X2,  . . . ,  Xn-i,  Xn),  via  the  formula 


9k{x) 


det  A(^Xi,  X2,  •  •  •  ,  Xk—2,  Xjg—i,  X,  Xjs-^-l,  Xk-^-2,  •  •  •  ,  Xn—l,  Xn) 

det  A(a;i,  X2,  ,  Xn-i,  Xn) 


(2) 


(here,  the  numerator  is  the  same  as  the  denominator,  but  with  x  in  place  of  Xk).  We  dehne 
D{xi,  X2,  . . . ,  Xn-i,  Xn)  to  be  the  modulus  of  the  determinant  of  A[xi,  X2,  . . . ,  Xn-i,  Xn), 
via  the  formula 


D{xi,  X2,  ...,  Xn-l,  Xn)  =  \  det  A{xi,  X2,  ...,  Xn-1,  Xn)  \.  (3) 

We  dehne  B  to  be  the  supremum  of  D{xi,  X2,  . . . ,  Xn-i,  Xn)  taken  over  all  sets  of  n  points 
Xi,  X2,  . . . ,  Xn-l,  Xn  in  S,  via  the  formula 

B=  sup  D{xi,  X2,  .  .  .  ,  Xn-l,  Xn).  (4) 

xi,X2,  ...,Xn-i,x„  in  S 

For  any  x  in  S,  we  dehne  u  =  u{x)  to  be  the  n  x  1  column  vector  dehned  via  the  formula 

Uk  =  fk{x)  (5) 

with  k  =  1,  2,  . . . ,  n  —  1,  n,  and  we  dehne  v  =  v{x)  to  be  the  n  x  1  column  vector  dehned 
via  the  formula 

Vk  =  9k{x)  (6) 

with  /c  =  1,  2,  . . . ,  n  —  1,  n. 

Theorem  2  below  asserts  the  existence  of  numerically  stable  n-point  interpolation  formu¬ 
lae  for  any  set  of  n  bounded  functions;  hrst  we  will  need  the  following  lemma. 

Lemma  1  Suppose  that  n  is  a  positive  integer,  S  is  an  arbitrary  set  containing  at  least  n 
points,  and  fi,  f2,  . . . ,  fn-i,  fn  a'k'G  complex-valued  functions  on  S  that  are  linearly  indepen¬ 
dent. 

Then,  there  exist  n  points  xi,  X2,  . . . ,  Xn-i,  Xn  in  S  such  that  the  vectors  u{xi),  u{x2), 
. . . ,  u{xn-i),  u{xn)  defined  in  (5)  are  linearly  independent. 

Proof.  We  apply  the  modihed  Gram-Schmidt  process  (see,  for  example,  [2])  to  the  set  of  all 
n  X  1  column  vectors  u{x)  dehned  in  (5)  for  all  x  in  S,  while  ensuring  that  all  pivot  vectors 

are  non-zero  via  appropriate  column-pivoting.  □ 
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Theorem  2  Suppose  that  S  is  an  arbitrary  set,  n  is  a  positive  integer,  fi,  f2,  ■  ■  ■ ,  fn-i,  fn 
are  bounded  complex-valued  functions  on  S,  and  e  is  a  positive  real  number  such  that 


e  <1. 

(7) 

Then,  there  exist  n  points  xi. 

X2,  . . . ,  Xn-i,  Xn  in  S  and  n  functions  gi,  g2,  ■  ■ 

•  ;  Qn—lj  Qn 

on  S  such  that 

\9kix)\  <  1  +  e 

(8) 

for  all  X  in  S  and  k  =  1,  2,  . . . , 

n  —  1,  n,  and 

n 

fix)  =  fixk)gkix) 

(9) 

k=l 


for  all  X  in  S  and  any  function  f  defined  on  S  via  the  formula 

n 

fix)  =  J2^kfk{x),  (10) 

k=l 

for  some  complex  numbers  ci,  C2,  ■  ■  ■ ,  Cn-i,  Cn- 


Proof.  We  assume  without  loss  of  generality  that  the  functions  /i,  /2,  . . . ,  /n-i,  fn  are 
linearly  independent. 

Then,  due  to  Lemma  1,  B  dehned  in  (4)  is  strictly  positive.  Since  the  functions  /i,  /2, 
. . . ,  fn-i,  fn  are  bounded,  D{xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn)  (dehned  in  (3))  is  also  bounded,  and 
hence  B  dehned  in  (4)  is  not  only  strictly  positive,  but  also  hnite.  Therefore,  there  exist  n 
points  Xi,  X2,  . . . ,  Xn-i,  Xn  in  S  such  that 


B  -  D{xi,  X2,  .  .  .  ,  Xn-l,  Xn)  <  —e  (11) 

and  D{xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn)  is  strictly  positive. 

Dehning  gi,  g2,  . . . ,  fi'n-i,  dn  via  (2),  we  obtain  (9)  from  the  Cramer  rule  applied  to  the 
linear  system 

Av  =  u,  (12) 

where  A  =  A[xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn)  is  dehned  in  (1),  v  =  v{x)  is  dehned  in  (6),  and  u  =  u{x) 
is  dehned  in  (5).  Due  to  the  combination  of  (11)  and  (7), 

-^  <  D{xi,  X2,  .  .  .  ,  Xn-l,  Xn),  (13) 


and  due  to  the  combination  of  (11)  and  (13), 


B 

D{xi,  X2,  Xn-1,  Xn) 


1  <c 


(14) 


we  also  observe  that,  due  to  (4), 


D{X\,  X2,  ■  ■  ■  )  Xk—2-)  Xk—lj  X,  Xfc_|_2,  ■  ■  ■  )  Xn—l,  Xn)  ^  B  (lb) 


for  all  X  in  S.  Now,  (8)  is  an  immediate  consequence  of  (2),  (3),  (15),  (14).  □ 
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Remark  3  Due  to  (8),  the  interpolation  formula  (9)  is  numerically  stable. 


Remark  4  When  calculations  are  performed  using  floating  point  arithmetic,  it  is  often 
desirable  to  “normalize”  the  set  of  functions  /i,  /2,  . . . ,  fn-i,  fn  before  applying  Theorem  2, 
by  replacing  this  set  with  the  set  of  functions  /i,  /2,  . . . ,  fn-i,  fn,  where  fk  is  the  function 
dehned  on  S  via  the  formula 


Mx) 


fk(x) 


(16) 


for  example. 


Remark  5  When  S  is  compact  and  the  functions  /i,  /2,  ...,  fn-i,  fn  are  continuous, 
Theorem  2  holds  with  e  =  0  rather  than  e  >  0,  since  a  continuous  function  D  on  a  compact 
space  S  attains  its  maximal  value.  For  the  same  reason.  Theorem  2  holds  with  e  =  0  when 
S'  =  R'^  for  some  positive  integer  d,  the  functions  /i,  /2,  . . . ,  fn-i,  fn  are  continuous,  and 
fk{x)  — >  0  as  |a;|  — >  cx)  for  all  A;  =  1,  2,  . . . ,  n  —  1,  n. 

The  following  theorem  formalizes  the  obvious  observation  that  integrating  both  sides 
of  (9)  yields  numerically  stable  quadrature  formulae. 

Theorem  6  Suppose  that  S  is  a  measure  space,  w  is  a  nonnegative  real-valued  integrable 
function  on  S  (that  serves  as  the  weight  for  integration),  n  is  a  positive  integer,  fi,  f2,  ■  ■  ■ , 
fn-i,  fn  are  bounded  complex-valued  integrable  functions  on  S,  and  e  <  1  is  a  positive  real 
number. 

Then,  there  exist  n  complex  numbers  wi,  W2,  ■  ■  ■ ,  vjn-i,  Wn  such  that 

\wk\  <  (1  +  e)  J  w{x)  dx  (17) 

for  all  k  =  1,  2,  . . . ,  n  —  1,  n,  and 

/n 

f{x)w{x)dx  =  ^Wkf{xk)  (18) 

k=i 

for  any  function  f  defined  on  S  via  (10),  where  xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn  are  the  n  points  in  S 
chosen  in  Theorem  2. 


Proof.  For  each  A;  =  1,  2,  . . . ,  n  —  1,  n,  we  dehne  Wk  via  the  formula 

=  J  9kix)  w{x)  dx,  (19) 

where  gi,  g2,  . . . ,  gn-i,  gn  are  dehned  in  (2).  Then,  (17)  is  an  immediate  consequence  of  (19) 
and  (8).  Moreover,  (18)  is  an  immediate  consequence  of  (9)  and  (19).  □ 


Remark  7  Needless  to  say,  the  weight  function  w  in  the  above  theorem  is  snperhuous;  it 
could  be  absorbed  into  the  measure  on  S.  However,  we  fonnd  the  formulations  of  Theorems  6 
and  10  involving  w  to  be  convenient  in  applications.  While  Theorems  6  and  10  require 
the  functions  /i,  /2,  . . . ,  fn-i,  fn  to  be  bonnded  (perhaps  after  “normalizing”  them  as  in 
Remark  4  or  otherwise  rescaling  them),  these  theorems  do  not  require  the  weight  function 
w  to  be  bonnded. 
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Remark  8  Theorem  6  asserts  the  existence  under  very  mild  conditions  of  numerically  sta¬ 
ble  quadratures  that  integrate  linear  combinations  of  n  functions  using  the  values  of  these 
linear  combinations  tabulated  at  n  appropriately  chosen  points.  In  contrast,  construction 
of  optimal  “generalized  Gaussian”  quadratures,  which  tabulate  the  linear  combinations  at 
fewer  nodes  than  the  number  of  functions,  requires  more  subtle  techniques  (see,  for  example, 
the  references  cited  in  [5]). 

Remark  9  The  proof  of  Theorem  2  does  not  specify  a  computational  means  for  choosing  the 
points  xi,  X2i  . . . ,  Xn-ii  Xn  SO  that  (11)  is  satished  (that  is,  so  that  the  resulting  interpolation 
and  quadrature  schemes  are  guaranteed  to  be  numerically  stable).  However,  combining  the 
algorithms  described  in  [1],  [3]  with  appropriate  discretizations  of  S  yields  methods  for 
choosing  the  points  that  are  proven  to  work,  both  in  theory  and  in  practice. 

Theorem  2  provides  the  bound  (8)  under  the  rather  weak  assumption  that  the  functions 
fi,  f 2,  ,  fn-i,  fn  are  bounded  (in  fact,  this  assumption  is  necessary  for  (8)).  The  following 
theorem  provides  a  stronger  bound  under  the  additional  assumption  that  there  exists  a 
measure  with  respect  to  which  the  functions  /i,  /2,  . . . ,  fn-i,  fn  are  orthonormal. 

Theorem  10  Suppose  that  n  is  a  positive  integer,  S  is  a  measure  space  containing  at  least 
n  points,  w  is  a  nonnegative  real-valued  integrable  function  on  S  ( that  serves  as  the  weight 
for  integration),  fi,  f2,  . . . ,  fn-i,  fn  ore  bounded  complex-valued  sguare-integrable  functions 
on  S,  and  e  <  1  is  a  positive  real  number.  Suppose  further  that  fi,  f2,  . . . ,  fn-i,  fn  are 
orthonormal,  that  is, 

j  \fk{x)\^  w{x)  dx  =  I  (20) 

for  all  k  =  1,  2,  . . . ,  n  —  1,  n,  and 

J  fj{x)  fk{x)w{x)  dx  =  0  (21) 

whenever  j  ^  k. 

Then, 


for  all  X  in  S  and  k  =  1,  2,  . . . ,  n  —  1,  n,  where  gi,  §2,  .  ■  ■ ,  Qn-i,  Qn  are  defined  in  (2),  with 
the  n  points  xi,  X2,  . . . ,  Xn-i,  Xn  in  S  chosen  in  Theorem  2. 

Proof.  In  order  to  prove  (22),  for  each  k  =  1,  2,  . . . ,  n  —  1,  n,  we  define  the  function  hk  on 
S  via  the  formula 

n 

hk{x)  =  Y.fji^) 

i=i 

and  demonstrate  that 


fj{y)9k{y)w{y)dy 


5 


and 


(25) 


hk{x)  =  gk{x) 

for  all  a;  in  S'  and  A;  =  1,  2,  . . . ,  n  —  1,  n. 

Suppose  that  /  is  defined  via  (10).  We  first  show  that 

n 

fi^)  =  (26) 

k=l 

for  all  X  in  S.  Substituting  (23)  into  the  right  hand  side  of  (26)  and  exchanging  the  orders 
of  summation  and  integration, 

n  n  p  n 

fi^k)  hk{x)  =  fj(x)  /  fj{y)  Y  fi^k)  gk{y)  w{y)  dy  (27) 

k=l  j=l  k=l 

for  all  X  in  S.  Due  to  the  combination  of  (27)  and  (9), 

Y  fi^k)  hk{x)  =  Y  fji^)  /  fjiv)  fiv)  ^(y)  dy  (28) 

k=i  j=i 

for  all  X  in  S.  Then,  (26)  is  an  immediate  consequence  of  applying  (10),  (20),  and  (21)  to 
the  right  hand  side  of  (28). 

We  now  complete  the  proof  of  (25).  For  any  x  in  S',  due  to  (26), 

At  =  u,  (29) 

where  A  =  A(xi,  X2,  ■  ■  ■ ,  Xn-i,  Xn)  is  dehned  in  (1),  m  =  u(x)  is  dehned  in  (5),  and  t  =  t(x) 
is  dehned  to  be  an  n  x  1  column  vector  via  the  formula 


tk  =  hk(x)  (30) 

with  A;  =  1,  2,  . . . ,  n  —  1,  n;  subtracting  (12)  from  (29), 

A(t-v)  =  t)^  (31) 

where  v  =  v(x)  is  dehned  in  (6).  Since  /i,  /2,  . . . ,  fn-i,  fn  are  orthonormal,  they  are  linearly 
independent.  Thus,  due  to  Lemma  1,  B  dehned  in  (4)  is  strictly  positive,  so  that  A  dehned 
in  (1)  is  invertible,  and  therefore,  due  to  (31), 


t(x)  =  v(x) 


(32) 


for  all  X  in  S'.  Then,  (25)  is  an  immediate  consequence  of  (32),  (30),  (6). 
Finally,  due  to  the  Cauchy-Schwarz  inequality. 


fk(y)9k(y)w(y)  dy 


<\  \fk(y)\^w(y)dy  J  /  \gk(y)\^  w(y)  dy, 


(33) 


and  due  to  (8), 

]J J  Idkiy)]"^  w(y)  dy  <  (1  +  e)  ^ J  w(y)  dy  (34) 

for  all  A;  =  1,  2,  ... ,  n  — 1,  n.  Then,  (24)  is  an  immediate  consequence  of  (23),  (33),  (20),  (34), 
and  then  (22)  is  an  immediate  consequence  of  (24)  and  (25).  □ 
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Remark  11  Due  to  (22),  the  interpolation  formula  (9)  is  numerically  stable.  While  the 
numerical  stability  guaranteed  by  (8)  is  sufficient  under  most  conditions,  sometimes  the 
bound  (22)  is  more  useful. 

Remark  12  Theorem  10  generalizes  easily  to  the  case  when  the  functions  fi,  f 2,  ,  fn-i,  fn 
are  not  precisely  orthonormal,  bnt  only  “close”  to  orthonormal,  in  the  sense  that  the  condi¬ 
tion  nnmber  of  their  Gram  matrix  is  reasonably  small. 

Remark  13  One  often  encounters  inhnite-dimensional  spaces  of  fnnctions  that  are  hnite- 
dimensional  to  any  specified  precision.  A  typical  sitnation  of  this  kind  involves  the  range  of 
a  compact  operator,  and  the  usual  way  to  construct  the  hnite-dimensional  approximation  is 
via  the  Singular  Value  Decomposition  (see,  for  example,  [5]).  When  combined  with  this  ob¬ 
servation,  the  apparatus  of  the  present  note  becomes  applicable  to  many  inhnite-dimensional 
spaces  of  functions. 
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